%This script estimates bounds of the identified set
%this replicates the code in distIVbounds.ado
%all the analysis is done for the appropriate bootstrap sample




%generate variables for appropriate subsample
D =regressors(sample,1);
Z =regressors(sample,2);
poor =regressors(sample,3);
schid =regressors(sample,4);
missing=regressors(sample,5);
missing=logical(missing);
notmissing=~missing;
X=regressors(sample,6:158);

%%
samplenotmissing=cumsum(notmissingoriginal);
samplenotmissing=samplenotmissing(sample(notmissing));

%%
IF1upper =IFbounds(samplenotmissing,1);
IF1lower =IFbounds(samplenotmissing,2);
IF0upper =IFbounds(samplenotmissing,3);
IF0lower =IFbounds(samplenotmissing,4);

%%
%basic calculation of bounds
beta= X(notmissing & Z==1,:)\IF1upper(Z(notmissing)==1);
ASF1upper= X *beta;
beta= X(notmissing & Z==1,:)\IF1lower(Z(notmissing)==1);
ASF1lower= X *beta;

beta= X(notmissing & Z==0,:)\IF0upper(Z(notmissing)==0);
ASF0upper= X *beta;
beta= X(notmissing & Z==0,:)\IF0lower(Z(notmissing)==0);
ASF0lower= X *beta;


gphiupper=max(ASF1upper-ASF0lower,0);
gphilower=max(ASF1lower-ASF0upper,0) ;

%%
bounds=[mean(gphiupper(poor==1)) mean(gphiupper(poor==0)) ...
    mean(gphilower(poor==1)) mean(gphilower(poor==0))];